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Abstract 

The subject of this study is the exchange-correlation-energy functional of reduced density matrix 
functional theory. Approximations of this functional are tested by applying them to the homoge- 
neous electron gas. We find that two approximations recently proposed by Gritsenko, Pernal, and 
Baerends, J. Chem. Phys., 122, 204102 (2005), yield considerably better correlation energies and 
momentum distributions than previously known functionals. We introduce modifications to these 
functionals which, by construction, reproduce the exact correlation energy of the homogeneous 
electron gas. 

PACS numbers: 71.10.-w,71.10.Ca,05.30.Fk,31.25.-v 
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I. INTRODUCTION 



Reduced-density-matrix-functional theory (RDMFT) is one possible way to tackle the 
problem of electronic correlation. It is based on Gilbert's theorem which is an extension of 
the Hohenberg-Kohn theorem to non-local external potentials. Gilbert's theorem guarantees 
that the ground-state expectation value of any observable of a quantum mechanical system 
is a unique functional of the ground-state one-body reduced-density matrix (1-RDM). Thus, 
the fundamental quantity in RDMFT is the 1-RDM instead of the electronic density on which 
DFT is built upon. Theproperties of the exact 1-RDM have been the subject of theoretical 
studies for a long time [l|, l2|, y, |j] . Nevertheless, only relatively recently approximate total- 



energy functionals of the 1-RDM were used in practical applications 
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15l . Il6l . Il7j . All these approximate functionals give a satisfactory account of electronic 
correlations in small atoms and molecules at the equilibrium distance. The latest generation 
of functionals performs equally well at the molecular dissociation limit. Especially this latest 
success makes RDMFT particularly appealing. 

A great advantage of RDMFT, compared to DFT, is that the exact many-body kinetic 
energy is easily expressed in terms of the 1-RDM. More specifically, the total energy of a 
many-electron system in its ground state, characterized by an external potential V^(r), can 
be expressed in terms of the ground-state 1-RDM, 7 as 

1 



^tot [7] = ^ion + / d'r / rfV5(r-r 



2 



7lr,r ^ 



+ j dh j dV(5(r-r')l^(r)7(r,r') 



+ 



d^r / d^r 



,7(r,r) 7(r',r') 



r — r 



+ E^c [7] • (1) 



The first term in the above expression is a constant not related to the electronic degrees of 
freedom, for example the ion-ion repulsion energy. The next three terms are, respectively, 
the kinetic, the external potential and the electrostatic energy and they are known, explicit 
functionals of 7. Finally, the last term is the exchange and correlation (xc) energy which 
contains all the remaining electronic contributions to the total energy. The exact form of 
this functional is unknown and for practical applications needs to be approximated. Most of 
the approximate functionals for the xc energy that have been introduced so far are implicit 
functionals of 7. They depend explicitly on the natural orbitals (pj, and the corresponding 
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occupation numbers rij which are defined as the eigenfunctions and eigenvalues of the 1- 
RDM: 

j rfV7(r,rO¥.,(r')=n,y,,(r). (2) 

Viewed as a functional of arbitrary 7 in an appropriately defined domain, the functional 
given in Eq. ([1]) has a minimum value at the ground-state 7. This appropriate domain 
is defined through subsidiary conditions for 7 known as A^-representability conditions. It 
was shown by Coleman p^] that there are two such conditions for 7, and they concern the 
occupation numbers 

00 

^ rij = , < < 1 , (3) 
i=i 

where N is the total number of electrons. These conditions guarantee that 7 corresponds 
to either a pure many-body state or an ensemble of pure states. The first condition can be 
enforced in the minimization with respect to 7 through the Lagrange-multipliers method. 
In that way, the quantity to be minimized becomes 

(4) 

where /i is the corresponding Lagrange multiplier. yU was shown to be equal to the chemical 
potential, i.e. the derivative of the total energy with respect to the total number of elec- 
trons 1^ . Interestingly, this allows one to exploit the discontinuity of as a function of the 
particle number for the calculation of the fundamental gap of materials and is a motivation 
for the development of RDMFT schemes for periodic systems in order to address questions 
like the semiconductor and insulator gaps. 

The second of the iV-represent ability conditions ([2]) has a dramatic consequence: it allows 
for border minima in the occupation number optimization. In other words, it allows for the 
the existence of a subset of the optimal occupation numbers which are exactly equal to either 
one or zero and do not satisfy the condition dF/duj = 0. We refer to the corresponding 
states as pinned states. It is rather unlikely for the exact theory to produce pinned states 
for most systems of interest. For an occupation number being exactly equal to one, the 
corresponding natural orbital would be present in all determinants of the full CI expansion 
with non-zero coefficients. A situation like that has not been found for small atoms and 
molecules or for the HEG, where the exact solution can be calculated. While this is true 
for the exact xc functional, all the approximate functionals which are considered in this 



work, yield pinned states with nj = 1 for all systems we applied them to, except for the 
two-electron systems. These pinned states are core states and, in the exact theory, they 
correspond to occupation numbers which are marginally smaller than one. Hence, as far as 
the optimal 7 is concerned, the approximate result, rij = 1, is perfectly satisfactory. The 
important implication, however, is that dF/ dnj 7^ for the pinned states, and consequently 
6F/6'y{r, r') 7^ at the optimal 7. 

A number of approximate functionals for E^c, including those of interest in the present 
work, can be cast into the form 



(5) 



(nj.ni) 
J] and 



i.e. they have the form of the usual Hartree-Fock exchange modified by the function / 
of the occupation numbers. The first such approximation was introduced by Miiller 
corresponds to the function f{nj,ni) = ^UjUi. Miiller considered a more general exponent 
for the occupation number product in the exchange-like term and found an optimal exponent 



of 1/2. By modelling the exchange and correlation hole, Buijse and Baerends [5j arrived at 
the same functional. Goedecker and Umrigar ^ (GU) considered a modification by explicitly 
removing the self-interaction (SI) terms. They also presented a direct minimization with 
respect to the natural orbitals and the occupation numbers and found correlation energies 
for small atomic systems which are in very good agreement with the exact results. Later 
however, it was realized that the GU functional fails to reproduce the correct dissociation 



limit for small molecules 
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2l| . On the other hand, the Miiller functional yields the correct 



dissociation limit but, in all cases, overestimates substantially the correlation energy 
the 
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ast decade, several other functionals of the 1-RDM have been introduced [8|, 19, 
3, 16 



Id, 



171] and applied to atomic and molecular systems. Recently, Gritsenko 



et al. [13(1 proposed improved 1-RDM functionals based on a hierarchy of repulsive corrections 
to the Miiller functional. In that way, they attempted to correct the overcorrelation of this 
functional. The functionals corresponding to these hierarchical corrections, are called BBCl, 
BBC2, and BBC3. For all these functionals, it is essential to divide the natural orbitals into 
strongly and weakly occupied ones. This distinction appears naturally for finite systems 
since usually a subset of the orbitals corresponds to occupation numbers close to one, and 
the rest to occupation numbers close to zero. For the BBCl and BBC2 functionals, the 



BBC2: f{n„ni) = { 



function f{nj,ni) is: 

„„„ „ , ; /.and ../weakly occupied, 

BBCl: f{nj,ni) = '{ (6) 

^Tij Til , otherwise, 

-y/rij rii , j 7^ /, and j, I weakly occupied, 

5 3 7^ h and j, / strongly occupied, (7) 
^rij rii , otherwise. 

Finally, in the BBC3 functional the anti-bonding orbital is treated as strongly occupied 
orbital. Additionally, the self-interaction terms are removed as in the GU functional, except 
for the bonding and anti-bonding orbitals. Gritsenko et al. 13| applied the BBC functional 
to diatomic molecules and showed that they give an accurate description of these molecules 
at both the equilibrium distance and the dissociation limit. 

There is a strong motivation for the extension and application of 1-RDM functionals to 
solid-state systems. This motivation stems from the success of these functionals in the de- 
scription of electron correlation for finite systems, as well as the difficulties of DFT methods 
in describing certain materials and properties such as the band gap of semiconductors and 



insulators 
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251 ] or the band width of the conduction band in Na 26|. 
A very important prototype system, which serves as a benchmark for the performance 
of approximate 1-RDM functionals when applied to periodic systems, is the homogeneous 
electron gas (HEG). Furthermore, this system can serve as a laboratory for the development 
of approximate functionals, in a fashion similar to DFT. As a consequence of translational 
invariance, the natural orbitals can be chosen as plane waves and the search for the ground 
state 1-RDM is restricted to the optimization of the momentum distribution n(k), which is 
the occupation number that corresponds to the plane-wave natural orbital with wavevector 
k. An important point to note is that the self-interaction terms for the plane-wave natural 
orbitals vanish. Consequently, the GU and the Miiller functionals are identical. Finally, 
approximations that involve a special treatment of single orbitals have zero effect in the 
continuous wavevector case. Thus, the BBC3 functional is identical to BBC2. 

As a consequence of the rotational invariance of the HEG the occupation numbers have 
the property n{k) = n{k), i.e. all the plane-wave natural orbitals corresponding to the 
same absolute value k are degenerate (with respect to the occupation number). This allows 
one to perform unitary transformations among the degenerate plane waves leading, e.g., to 



angular-momentum eigenfunctions 

Vkim{r) = ji{kr) Yirai^) (8) 

for the natural orbitals, where ji{kr) are spherical Bessel functions and Yim spherical har- 
monics. Since these functions are localized in real space it is conceivable to include self- 
interaction corrections in terms of the natural orbitals ([8]) for the HEG. To our knowledge, 
this possibility has not been explored so far and is also beyond the aim of the present work. 
Choosing the natural orbitals as plane waves, the 1-RDM of the HEG can be written as 

k 

where V is the volume of the system (V oo). Substituting this expression in Eqs. 
and ([5]) we obtain for the total energy 

^tot = 2 5^ f n{k,) -^J2f (^(^i) ' ^(^2)) , (10) 

1^ ^ ^ |ki-k2| +«2 

where ki, k2 are wavevector indices, and a is a small quantity (usually ~ 10~^) included 
for numerical stability. As in Hartree-Fock, the external potential energy and the electronic 
Coulomb repulsion energy, i.e. the third and fourth terms in Eq. ([1]) cancel exactly with the 
first term i.e. the ion-ion interaction. The quantity F of Eq. (j4]) per particle then becomes 



F 3 



N 2kl JO 



dki kl {kl - 2fi) n{ki) 

{kl - k2f + 



+ ,_,^3 / dki I dk2 kl k2 log 



4vr/i;^ JO JO 



_{ki + k2f + a^ 

where kp = (97r/4)^/^r~^ is the Fermi- wavevector of the non- interacting HEG and is the 
radius (in atomic units) of the sphere with volume equal to the volume per electron. 



f{n{ki),n{k2))+fi, (11) 



Cioslowski and Pernal 



27| applied the Miiller functional to the HEG and calculated 



analytically the resulting momentum distribution 

n{k) = 51271 p{l + 4k'^)-\ (12) 

where p is the electron-density-per-spin, p = 3(87rrg)~^. The corresponding total energy 
per particle is independent of the density and equal to —1/8 Hartree. It is obvious that the 
solution of the Eq. (fT2|) is consistent with the second iV-represent ability constraint of Eq. ([3]) 



only for p < (5127r)~^, i.e. > 5.77. In other words, the Miiller functional gives a solution 
without pinned states only for > 5.77. In addition, Cioslowski and Pernal demonstrated 
that the Oxford-Lieb 28| bound is violated for p > 1.65 x 10^^, i.e. < 4.167. The solution 



with pinned states for < 5.77, was calculated by Csanyi and Arias [7]. More specifically, 
for rg < 5.77, one gets an optimal momentum distribution n{k) with n{k) = 1 for k below a 
certain value kp and fractional n{k) < 1 for A; > kp. This behavior is in complete analogy to 
the case of finite systems for the Miiller functional. Unfortunately, it is in conflict with the 
fact that the exact momentum distribution 29|, ISOj is a monotonically decreasing function 



of k and is strictly smaller than 1, i.e. there are no pinned states. Additionally, the exact 
momentum distribution is concave for k < kp, it shows a discontinuity at kp, and for k > kp 
it goes to zero asymptotically. The size of the discontinuity is decreasing with r^. 

In addition to the Miiller functional, Csanyi and Arias [7] considered a similar functional 
derived from a tensor product expansion of the two-body density matrix, which they called 
Corrected Hartree-Fock (CHF). Unfortunately, CHF gives zero correlation for the HEG in 
the high-density limit (r^ — > 0), coinciding with Hartree-Fock. In the opposite limit, it 
strongly overcorrelates giving the same results as the Miiller functional. In the intermediate 
region, including the metallic densities, the result for the correlation energy is close to 
the exact but its dependence on is monotonically decreasing instead of increasing. In 
an attempt to improve over the Miiller functional and CHF, Csanyi, Goedecker and Arias 



considered an improved tensor product expansion of the two-particle density matrix [8|. The 
resulting functional, which is called CCA, performs very well in the high-density regime and 
significantly better than the previous two functionals in the region of metallic densities. The 
deviation from the exact correlation energy increases with and at higher densities CGA 
coincides with CHF and the Miiller functional. 

In the present work, we apply the BBCl and BBC2 functionals of Gritsenko et al. [l^ 
to the HEG and compare with previous functionals as far as the resulting correlation ener- 
gies are concerned. We also investigate other features like the resemblance of the resulting 
momentum distribution to the exact and the state-pinning. In order to apply the BBCl, 2 
functionals to the case of the HEG through Eq. ffTTl) . we need to distinguish between strongly 
and weakly occupied orbitals. For finite systems, Gritsenko et al. 13| chose the first A^/2 
natural orbitals to be strongly occupied. In complete analogy, we can use a critical wavevec- 
tor kc = kp below which all states are assumed to be strongly occupied while above they are 
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weakly occupied. 

An additional goal of the present work is to demonstrate that the HEG can be used to 
develop functionals suitable for metallic systems. The idea is to modify the BBCl functional 
in such a way that it yields the exact correlation energy for the HEG. This is achieved in two 
different ways: (i) For each given density we choose such that BBCl reproduces the exact 
correlation energy of the HEG at that density. We call this functional fcc-functional. (ii) 
We introduce a function s{rs) multiplying the xc terms of Eq. ([5]) for two weakly occupied 
orbitals, keeping kc = kp. In this way, we replace the sign change of the BBCl functional 
with the parameter s. Accordingly, we call this functional s-functional. 

In the following section we present details of the numerical implementation as well as the 
results of applying the BBCs and the kc- and s-functional to the HEG. 

II. NUMERICAL IMPLEMENTATION, RESULTS 

The minimization of the energy expression ( ITT]) with respect to n{k) is performed using 
the steepest-descent method. We choose to work in energy-space instead of fc-space, i.e. we 
perform the variable substitution e = k'^/2 and solve numerically the minimization problem 
for n{e). Working in energy-space rather than /c-space improves the stability of the numerical 
treatment. The energy e is discretized using a double- logarithmic mesh centered at ep = 
kp/2, where the occupation varies the most. The upper limit of integration is chosen such 
that the momentum-distribution function has dropped to values smaller than 10~^. The 
double integration with respect to the energy is carried out using an adaptive grid technique 
capable of treating integrable singularities like the logarithmic singularity of the present 
problem. The values of the energy-distribution function n{e) in between the mesh points, 
necessary for the adaptive grid method, are obtained from an interpolation scheme. The 
A^-representability constraint < n{e) < 1 is implemented through the substitution n(e) = 
sin^[7r 9{e)/2] and variation with respect to 0{e). Extra care is required to avoid n{e) being 
falsely pinned to or 1. Indeed, if for a particular e, n(e) gets very close to 1 or during the 
variation it would stay pinned at that point. The variation with respect to the Lagrange- 
multiplier fi is implemented as an external iterative procedure, thus achieving convergence 
for each value of fi. The correct value of fi is selected by requiring the momentum-distribution 
function to integrate to the correct number of electrons for a given value of Vg. Finally, we 
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FIG. 1: (Color online) The correlation energy of the HEG as a function of calculated with 
the BBCl and BBC2 functionals compared to various other calculations. The Monte Carlo result 
represents the Perdew-Wang fit 



3ll | of the DMC data of Ortiz and Ballone 



corresponds to the numerical results by Csanyi and Arias 



291]. The dotted line 



71] employing the Miiller functional, for 



Vs < 5.77. Its continuation, the dash-dotted line, for > 5.77, stands for the analytical results of 



Cioslowski and Perna! 



as well as the CGA 



271 ] employing the same functional. The results for the CHF functional 



3] are also shown. 



found that a reasonable value for the parameter a in Eqs. f llUI) and fill I) is 10 ^, which avoids 
both numerical problems as well as the dependence of the results on a. 



A. Application of the BBC functionals to the HEG 

In Fig. [H we show the correlation energy of the HEG as a function of r^. The correlation 
energy calculated with the BBCl and BBC2 functionals is significantly closer to the exact 
than any other functional over the whole range of r^. Both functionals also seem to reproduce 
the correct asymptotic limit of zero correlation for the dilute HEG. For small densities up 
to metallic densities, the BBC functionals under-correlate, i.e. the absolute value of the 
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correlation energy is too small. In the dense limit, they over-correlate and the crossover is 
at around rs=0.5 and 0.3 for the BBCl and BBC2, respectively. Unfortunately, in the area 
of low metallic densities both functionals yield correlation energies which deviate from the 
exact values by 50%. In absolute numbers, the error of BBCl and BBC2 is of the same 
order as the RPA result 3l|]. Nevertheless, in the range 0.1 < < 1, the BBC functionals 



perform remarkably well. Compared to all previous 1-RDM functionals, BBCl and BBC2 
offer a much better account of the electron correlation for the HEG. Although less accurate 
than the CCA in the high density region, they perform better for metallic densities and they 
reproduce the limit of zero correlation at the dilute HEG limit where the MiiUer functional, 
CHF and CCA fail. 

A feature of the exact momentum distribution, namely the discontinuity at the Fermi 
wavevector kp, is reproduced by the BBC functionals. The discontinuity is more pronounced 
for the BBCl functional, as can be seen in Fig. [21 where we plotted the momentum distri- 
bution of the HEG with = 1 and = 5 using the Miiller functional as well as BBCl and 
BBC2. Contrary to BBCl and BBC2 the Miiller functional does not yield a discontinuity. 
To our knowledge, there is no report of any other 1-RDM functional reproducing this fea- 
ture of the HEG. To extract the size of the discontinuity quantitatively, we used two energy 
mesh points very close to ep (at a distance of ±10~'^eF). In Fig. [3l we plot the size of the 
discontinuity An as a function of r^. As we see, it increases monotonically with and has 
the tendency to saturate for large for both BBCl and BBC2. For BBC2 the discontinuity 
is substantially smaller than for BBCl and it goes to zero at 0.6 < < 0.7. Both the size 
and the dependence on are in complete disagreement with the exact theory, where An 
is substantially bigger and decreases with as one can see from the two fits to the DMC 



data 



29 



30|. 



We now turn to the question of state pinning. As we see in Fig. [2l state pinning is 
a common feature of all the functionals we employed. To verify that the states are truly 
pinned, we plot the functional derivative 6F/6n{k) which is non-zero for pinned states. 
As one expects, the number of pinned states decreases with increasing r^. We define a 
wavevector kp below which the corresponding states are pinned, i.e. n{k) = 1 for A; < kp. 
In Fig. m we plot kp as a function of for the Miiller and the two BBC functionals. 
For all three functionals, kp decreases monotonically with r^. For the Miiller functional our 
numerical calculation confirms the analytic result 271] that above a critical value of = 5.77 
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FIG. 2: (Color online) The momentum distribution n{k) of the HEG for rs=l (left) and = 
5 (right) calculated with the Miiller functional, BBCl, and BBC2. BBCl and BBC2 show a 
discontinuity of the momentum distribution at kp. For comparison we include the fit to the DMC 
data of Ortiz-Ballone [29}]. The derivative 6{F/N)/6n{k) is also plotted. Note that for small k the 
derivative is not zero and n{k) is pinned at one. 



there are no pinned states, and therefore kp goes to zero at this value. This can already be 
seen in Fig. |2l where, for the Miiller functional at = 5, the derivative is very close to zero 
even for small wavevectors. Interestingly, for BBCl and BBC2 we found no such critical 
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FIG. 3: (Color online) The discontinuity An of the momentum distribution at fep for BBCl and 
BBC2 as a function of r^, compared to the fits of Ortiz and Ballone (OB) [2^, and Gori-Giorgi 



and Ziesche (GZ) [30( to the DMC data. 



value up to = 20. Indeed, for both BBCl and BBC2 the decrease of kp is much smaller 
than for the Miiller functional, kp being almost constant for BBC2. 

As we have seen, the performance of the BBC 1,2 functionals for the HEG is improved 
significantly compared to previous functionals. This is especially remarkable, given that 
they were originally constructed to describe the dissociation of small molecules. 



B. Improved functionals for the HEG 

In this section, we attempt to improve over the BBC functionals for the HEG. We present 
two functionals which are simple modifications of the BBCl functional. In both cases, the 
BBCl functional is generalized by introducing an additional function of such that the 
correlation energy of the HEG is reproduced exactly for each r^. As exact results we regard 
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FIG. 4: (Color online) The wavevector kp, where the optimal momentum distribution becomes 
fractional, as a function of r^. 



the Perdew-Wang fit 31i of the correlation energy obtained from diffusion Monte Carlo 



29|. The two Monte 



(DMC) calculations by Ceperley and Alder [3^ and Ortiz and Ballone 
Carlo calculations yield almost identical correlation energies. 

For the first functional, we adjust the critical wavevector kc, which is used to distinguish 
between strongly and weakly occupies states, instead of as in BBCl. We call 

this functional the /cc-functional. The corresponding function / in Eq. ([5]) then reads 



f{n{k,),n{k2)) 



^Jn{ki) n(k2) , otherwise. 



(13) 



We perform the fitting of kc over the range of metallic densities, 0.5 < < 5. The results 
are compiled in Table [H For = 0.5, kc ~ ^f, since this point is almost exactly reproduced 
by BBCl (see Fig. [T]). For > 0.5, kc is a monotonically increasing function of r^. Fitting 
kc has a strong impact on the momentum distribution which is displayed in Fig. O for 
Ts = 1 and Ts = 5. It is not surprising that the discontinuity is displaced from kp to kc- 
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TABLE I: The the critical wavevector kr for some metalhc densities. 
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FIG. 5: (Color online) The momentum distribution for the fcc-functional compared with the BBCl 
and the fit to the DMC data of Ortiz-Ballone |29j], for = 1 and = 5. For the /cc-functional the 
discontinuity is moved from kp to kc, i.e to 1.032A:f and 1.172A:f respectively. 

Additionally, its size is reduced significantly compared to BBCl. Both the displacement of 
the discontinuity and the decrease in the step size are in disagreement with the exact result. 

An alternative idea is to keep = kp fixed and consider a fitting parameter s multiplying 
the exchange-like terms when both states, ki and k2, are weakly occupied, i.e. 



finik,),nik2)) 



(14) 



-s{rs) \/n{ki) n{k2) , ki,k2>kp, 
\Jn{ki) n{k2) , otherwise. 

In this way, the parameter s which is a kind of strength of the xc terms becomes a function 
of as indicated in Eq. f ll4p . We call this functional the s-functional. The values of s 
for the fitting to the Ortiz and Ballone DMC results are included in Table [TTl As we see, 
s varies between 4 and —0.26 for the range of densities, O.l < < 10, we considered. In 
Fig. [6], we show the resulting, monotonically decreasing function s{rs). 

From Fig. [7] one can see that one of the advantages of fitting s instead of k^ is that the 
discontinuity of the momentum distribution remains fixed at kp. In addition, its size is 
almost constant (~ 0.2) as a function of r^. As we see in Fig. |3l the exact discontinuity is 
significantly higher and it is a decreasing function of Tg. Therefore, concerning the size of 
the discontinuity, the s-functional does not improve over the BBCl, which, in the dilute gas 
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0.1 


4.913 


1.0 


0.435 


0.2 


2.751 


1.5 


0.190 


0.3 


1.867 


2.0 


0.059 


0.4 


1.390 


3.0 


-0.074 


0.5 


1.087 


4.0 


-0.146 


0.6 


0.877 


5.0 


-0.189 


0.7 


0.727 


7.0 


-0.234 


0.8 


0.602 


10.0 


-0.263 



TABLE II: The fitted values of s for various values of for the s-functional. s was fitted to 
reproduce the correlation energies of the DMC calculation of Ortiz and Ballone 



limit, is close to the exact result. However, the increasing behavior of BBCl is improved by 
the s-functional which yields discontinuities that are almost constant as a function of r^. 

In conclusion, the s-functional results in momentum distributions that resemble the exact 
ones over the whole range of more closely than any of the other functionals considered 
here. 



III. CONCLUSION 

We have applied a variety of 1-RDM functionals to the HEG. We show that the BBC func- 
tionals [l^ yield a significant improvement over previous functionals as far as the correlation 
energy is concerned. In addition, they yield a discontinuity of the momentum distribution 
at the Fermi wavevector in resemblance of the exact HEG theory. However, the size and 
the dependence on the density of this discontinuity are not in agreement with the quantum 
Monte Carlo results. 

By introducing an appropriately fitted function of Vg in the BBCl functional, we demon- 
strate that the exact correlation energy of the HEG can be reproduced with a smooth and 
monotonic fitting function. For this function, we either use the critical wavevector /cc(^s) 
which distinguishes between the strongly and weakly occupied states, or a strength s{rs) 
multiplying the exchange-like terms for two weakly occupied states. Both of these proce- 
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FIG. 6: (Color online) The dependence of the fitting parameter s on resulting from fits to 



the Perdew-Wang parameterization 
results: the Ceperley and Alder 



3l[ | of the correlation ener gy f rom two different sets of DMC 



32(1 and the Ortiz and Ballone 



1 

0.8 
-0.6 
0.4 
0.2 
0, 



' 1 ' 


! ' 

^ r=l - 






1 " 




j — ^-functional _ 




j - - BBCl 




! — DMC 


1 


L . - 







0.5 



1 



1.5 



1 

0.8 
0.6 
0.4 
0.2 
0, 



' 1 ' 






r=5 - 


^ ^ \ I 


s 


\ 1 


^-functional _ 


V. 


- BBCl 


y 


-- DMC 






1 


■■~-:rr:r 1_ 







0.5 



1 



1.5 



FIG. 7: (Color online) The momentum distribution given by the s-functional, compared with the 



29!], for rs 



1 (s 



0.435), and 



BBCl and the fit to the DMC data of Ortiz-Ballone 
(s = —0.189). The discontinuity remains at kp and is approximately equal to 0.2. 

dures were applied to the BBCl functional. The two functionals, yielding by construction 
the exact correlation energy of the HEG, are assessed by the quality of the resulting mo- 
mentum distributions. We show that choosing the second procedure, i.e. the s-functional, 
is superior to fitting kc- The discontinuity An of the momentum distribution resulting from 
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the s-functional is nearly constant as a function of and hence represents a significant 
improvement over BBCl and BBC2. However, the momentum distribution obtained by the 
s-functional still deviates significantly from the exact one. To remedy this, more complicated 
strategies have to be considered, possibly with the introduction of more fitting parameters. 

Our functional of choice, being derived from the HEG, is expected to yield good results 
for metallic systems. The application to finite as well as non-metallic periodic systems is 
not straightforward because is not well-defined in these cases. Hence, the necessity arises 
to map Ts onto other quantities characterizing these systems, or to involve a LDA-type 
prescription relating s{rs) to the local density. 
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